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In this paper we describe the rescattering process in optical field ionization through a one- 
dimensional model, which improves the well-known quasistatic model by adding the smoothed 
Coulomb potential in its second step. The above-threshold ionization spectra and high-order 
harmonic generation are calculated from this model. They are qualitatively in agreement with 
the quantum results and experiments. In particular, we find that this model is characterized 
by chaotic scattering. Modern nonlinear theory is used to analyze this dynamical system. 
It is found that singular self-similar fractal structure exists in the phase-dependence energy 
spectra, and the unstable manifolds constitute the chaotic scattering pattern. Our results 
also demonstrate a close connection of irregular trajectories with the energy spectra and 
high-order harmonic generation. We conclude that the chaotic behavior plays an important 
role in this one- dimensional model of optical field ionization. 
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I. INTRODUCTION 



In recent years there have been significant advances in high-intensity laser technology. 
One major use of these laser systems has been for studies of the response of atoms and 
molecules to such intense fields. It has led to the discovery of a whole range of nonper- 
turbative phenomena, such as, multiphoton ionization (MTI), above-threshold ionization 
(ATI), and high-order harmonic generation (HOHG). Because the fields are very intense, 
the traditional perturbation expansion of the wave function in terms of the field-free states 
will fail completely to describe the dynamics of the system. Therefore, these problems are 
challenging theoretically and have recently received much attention both experimentally and 
theoretically J 1 ' 2 ' 

It has been concluded by numerical calculation of the one-electron Schrodinger equa- 
tion that these phenomena can be understood in the context of single-electron ionization 
dynamics. Various nonperturbative theories are developed to explain ATI and HOHG phe- 
nomena. Among these theories the quasi-static model is commonly used in treatment of 
the optical field ionization (OFI) problem in the regime of low field frequency.' 3,4 ' In this 
model, an electron can be considered to acquire its energy in a two-step process: first the 
electron is removed from the atom, thereby overcoming the ionization potential, and then 
it interacts with the laser field. The ionization probability of the atom in the first step is 
usually described by the tunneling theory. Since it is in the continuum state, the evolution 
of an electron wave packet after tunneling is usually described by classical motions governed 
only by the laser field. This simple two-step quasi-classical model accounts for many of the 
detailed characteristics of the photon and electron emission process and provides significant 
insights into the excitation dynamics under intense field conditions. 

However, in the case of an intense linearly polarized laser, the amplitude of the quiver 
motion of an electron after tunneling is large. Then the probability that an ionized electron 
returns to the vicinity of the nucleus and is rescattered by the ion is not negligible. This 
problem attracts much attention and becomes an active point recently. Many classical and 
quantum theories are developed to treat this rescattering process J 5,6 ' 7 ' It is found that this 
process may increase the fraction of the electron with higher energy , and is responsible for 
other important phenomena , such as the 'cut off' law in high-order harmonic generation. 
Different from the above works, in this paper we try to consider the rescattering process 
in the frame of one-dimensional (ID) quasiclassical model, which is derived from natural 
extending of the original quasistatic model by adding a smoothed Coulomb potential in 
its acceleration step . In sec. II we will introduce this ID model. According to this model 
the above-threshold energy spectra and high-order harmonic generation in the rescattering 
process are calculated in sec.V. Our results show that the rescattering process increases 
the higher energy electron and generates high-order harmonics. Additionally, the 'cut off' 
law is identified by this model. In particular, based on nonlinear theories, we discover 
a new physical process - chaotic scattering underlying optical field ionization and give a 
complete dynamical explanation for them (sees. Ill and IV). We find that self-similar and 
fractal structure exists in the plots of ATI energy vs initial field phase. We also identify 
those unresolved regions of ATI energy with chaotic layers. The location and width of those 
singular regions are estimated by a nonlinear theory. In sec. VI we analyze the dynamical 
structure in Poincare section which causes the chaotic scattering. Finally in sec. VII we 
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summarize our main results and discuss some problems for future study. 



II. THE MODEL 

First let us recall the well-known quasistatic theory briefly. As 7 = (Io/21/p) 1 ^ 2 <C 
l,huj/2U p <C 1 and field strength F < F th = Z 3 jVon\^a.u. ( atomic units are used in 
this paper), tunneling ionization occurs. Here, Jo is the ionization potential of the electron, 
Up = f^f-ff is the ponderomotive potential, F t h is the threshold field, Z is the core charge 
and n e ff is the effective principle quantum number. The corresponding threshold intensity 
of F th is 

I th = 1.37 x 10 14 (/o/13.6e\/) 4 /Z 2 (W/cm 2 ). (1) 

Based on the Landau tunneling ionization theory, ' 8 1, Ammosov, et al.' 9 ' suggested a simple 
formula for the ionization rate of complex atoms: 

r(f) = 1.61^( 1 ^) 2 -"-'- 5 exp(-|/ ? ) a.u. (2) 

This formula is commonly used in estimating the ATI energy of the ionized electrons J 10 ^ 

The second part of the quasistatic procedure uses classical mechanics to describe the 
evolution of an electron packet. After tunneling, the electron motion in the field is given by 

x = x [cos(ujt)] + v 0x t + x 0x , 

v x = -v sin(ujt) + v 0x , (3) 

where v = ZeF/m e u,x = ZeF/m e u 2 . The energy associated with the velocity v 0x consti- 
tutes the ATI energy and can be evaluated by averaging over the external field 

£4 T/ = f/ p [l + 2sin 2 o ], (4) 

where 0o = denotes the field phase. 

In a linearly polarized intense laser field , the amplitude of the quiver motion is large. 
Then the ionized electrons are driven by the laser field to return to the vicinity of the 
ions and are rescattered by them. This process will result in many important phenomena, 
such as 'plateau' in ATI peaks, 'ring' in angular distribution and 'cut-off' law of high-order 
harmonic generation. To take the rescattering process into consideration, we improve the 
above model by adding a smoothed Coulomb potential in its second procedure. Then the 
governing equation is (the atomic units is used for simplicity), 

H = V -- -T== + exFcos(cut), (5) 

where a — 0.1 the 'softening' factor to normalize the Coulomb potential. 

The initial position of an electron borne outside the potential barrier at time to is deter- 
mined by the tunneling ionization theory, 



3 



U eff = 1 + eFx cos(ujt ) = I , 



(6) 



where Jo = 0.5 is the ionization potential. 

Take n/2 < ut < 3n/2 as an example. It is observed that the electron corresponding 
to the maximum field strength (ujt = it) is closest to the nucleus. As Lut — | or ^, the 
initial electron position tends to positive infinity. Due to the fact that an electron will be in 
the continuum state after tunneling ionization, we prefer to choose the initial energy rather 
than its initial momentum of an ionized electron born at time to to be zero. That is, if 
the external field is removed away at time to the ionized electron will escape away from the 
Coulomb potential. On the contrary, the condition that the initial momentum is zero cannot 
guarantee the escaping of the ionized electron J 11 ' Then the initial momentum of an ionized 
electron can be derived from the above assumption. These initial conditions are sufficient 
for us to study the rescattering effects in the following sections. 



III. SELF-SIMILAR FRACTAL STRUCTURE AND SCALING LAW 

It is convenient to introduce the compensated energy E c advocated by Leopold and 
PercivalJ 12 ' 



When an electron is ionized completely and cannot return to the nucleus again, the Coulomb 
potential is weak enough and E c tends to a positive constant value which is just the ATI 
energy (Eati) in ultrashort pulse laser. 




(7) 
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FIGURES 
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FIG. 1. Plots of electron energy vs field phase. Dashed lines and solid 

lines denote the results from original model and the present model respectively. a) 
F = 0.06a.u.,LO = 0.07a.u. b)F = 0Ma.u.,uj = 0.04242a.u.. 
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FIG. 2. Trapping time vs field phase. F = 0.06a.u.,u; = 0.04242a.u.. 

In Fig.l we demonstrate the phase dependence of the ATI energy. A Runge-Kutta al- 
gorithm of the 4th order in x and p is employed in the calculations. One thousand points 
distributed equally in the range [7r/2,37r/2] are chosen as the initial phase. It was observed 
that the curves calculated from the original model where the Coulomb potential is neglected 
are smooth over all the range; however, in our new model things are much different. The de- 
pendence of E ATI on the initial phases is poorly resolved in the region [3.32, 4.05], [4.52, 4.71] 
for u = 0.07a.u. and in the region [3.30,3.70], [4.69,4.71] for u = 0.04242a.w. Magnification 
of the unresolved regions shows that an arbitrarily small change in the initial phase may 
result in a large change in the final electron energy We define the escaping time of an 
electron as the moment T ^ satisfying 
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E c (t = T) = 0, E c (t > T) > 0. 



(8) 



Naturally, the trapping time of an electron is defined as AT = T—t . A comparison between 
Fig.l and Fig. 2 shows that the unresolved regions in Fig.l coincide with the singularities of 
the trapping time. 

We call a value S a singularity if, in any small neighborhood domain there is a pair of 
which have different signs in their final moments. These singular values relate to the infinite 
trapping time . Generally in phase space they constitute a nonattracting hyperbolic invariant 
set. This invariant set has a fractal structure. Figure 3 shows successive magnification of 
the unresolved region and demonstrates clearly the existence of self-similarity. 
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FIG. 3. Successive magnification of the plot for F = 0.06a.u.,uj = 0.07a.u.. 



To calculate the fractal dimension of the singular points (j> 8 , we can employ the uncer- 
tainty exponent technique to obtain the fractal dimension of the singular set J 14 ' We randomly 
choose many values of in an interval containing the fractal set. We then perturb each value 
by an amount e and determine whether the final momentum corresponding to initial 0, — e 
and + e have the same sign. If so , we say that the value is e-certain; if not, we say it 
is e-uncertain. We do this for several e and plot on a log-log scale the fraction of uncertain 
values /(e). The result is plotted in Fig.4 which shows a good straight line and indicates 
a power law dependence /(e) ~ e 7 , where 7 = 0.15. The exponent 7 is related to the 
dimension of the fractal set of the singular values by 

D = 1 - 7 = 0.85. (9) 

In another aspect, it is interesting to make an intensive investigation on the escaping 
process of an electron after tunneling. We find that the trajectories demonstrate much 
different behavior for different initial phase. Some electrons drifts away from nucleus without 
a collision, others can undergo multiple collisions. In Fig. 5 we plot the probability vs. the 
collision times by making a statistics on 4000 trajectories corresponding to different phases. 
It is found that, approximately the probability decays with increasing collision times and 
satisfies the following scaling law 
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Pool ~ e , (10) 

where k = 0.21 , N denotes the collision number. 

Obviously, the possibility of multiple collisions decreases exponentially with increasing 
collision number. Although chaos implies the existence of infinite times of collisions, its 
probability is zero. This is because the hyperbolic invariant set corresponding to the trapped 
trajectories has zero measure. 




FIG. 4. Calculation of fractal dimension in the case F = 0.06a.u.,u; = 0. 07a.M. 
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FIG. 5. Probability of multiple collisions obtained from making statistics on 4000 trajectories 
corresponding to different initial phase (circle) . The straight line denotes the scaling law mentioned 
in the text. 



IV. UNRESOLVED REGIONS RESULTING FROM CHAOTIC LAYERS 

In previous sections we discovered that the unresolved regions have an irregular, self- 
similar fractal structure. In this section we will discuss the dynamical source of these phe- 



7 



nomena. We will determine the location and width of the unresolved regions by using 
modern nonlinear theory. 

First let us review some important conclusions ' 15 ' on a two-dimensional (2D) area- 
preserving mapping described by 

E = E + AE, <j) = (j) + f(E) (11) 

In the neighborhood of a separatrix (E s ), where the period f(E s ) tends to infinity, even 
small changes in energy E may result in considerable changes in phase 0. This is the very 
cause of local stochasticity. The following expression can be served as a good evaluation of 
the stochastic region boundary 

lg-n>„. , 12 , 

In the rest of this section we shall employ the Kramers- Henenberger (or acceleration) gauge 
and atomic units will be used to simplify our formula. The variable transformations are 

P = p + —sinujt, X = x %-cosut. (13) 



The new Hamiltonian is 



2 J(X + ^cosujtf + a 



The corresponding equations of motion are 

• X + -^coscut , 

X = P, P = ^ r . 15 

(X + ^cosut) 2 + a)i 

The fundamental property of the Hamiltonian system is that the Cartan differential 2-form 

dP A dX — dH A dt (16) 
is preserved under Hamltonian flow. As the term -^cosut in (16) vanishes , the energy 

E ^ P ^-7m=« (17) 

is an integral constant. Otherwise, its evolution is time-dependent and can be described 
approximately by virtue of Eqs. (16) and (18) 



^ = —PXcosujt/{X 2 + af. (18) 
dt to 2 



To evaluate the jump in energy, we consider the following straight-line trajectory as the 
zerot/i-order approximation, 

P = v c , X = v c (t-r), (19) 
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where v c and r are the momentum and time, respectively , at X = 0. Then the jump of the 
energy for each return is obtained by integrating Eq.(19) from negative infinity to positive 
infinity, 

2g f+OO 

AE = — v 2 c (t - T)cosuot/(vl{t - r) 2 + afdt. (20) 

CO J— oo 

The above integral can be rewritten in a simpler form, 

. _, e f+°° cos(s) , 

AE = —smcor / -r— ^i-^ds, (21) 

l£ J-oo S z + Ot(jJ z /Vc 

The integral included in the above formula can be expressed in term of a Bessel function, 

^ f °° C ° St tft (22) 



Ki(z) 



where z = a^Lojv c Now we divide our problem into the following two cases. 
A. Case 1, E sl = 0. The jump of energy is rewritten as 

AE = 5E sl sin(f), (f> = cut, (23) 

where 

«„4r rrl> (24) 

and f c i is derived from Eq.(18) by setting E = E sl and X = 0. As the 'softening' factor a 
is usually chosen as a small value, the Keplerian Law of the two-body problem is available 
approximately. Therefore, we obtain the corresponding change of phase as follows, 

= 0+-^r- (25) 
l^ol 2 

Equations (24) and (26) construct a 2D mapping which reflects the dynamical properties of 
the system (15) in the vicinity of E — E sl . The area-preserving nature of this mapping, a 
consequence of the invariance of the Cartan differential 2- form (17), can be checked easily. 

In terms of the nonlinear theory mentioned at the beginning of this section, we conclude 
that there exists a chaotic layer in the vicinity of E = E s i, the width of it is determined by 
condition (13). Then we get 

0<E sl -E <(^io5E sl )l (26) 
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Initial Phase (rad) 



FIG. 6. The chaotic layers predicted theoretically, showing coincidence with those unresolved 
regions. F = 0.06cl.u.,uj = 0.07a.u. 

In Fig. 6 we plot the initial energy vs the initial field phase. The above inequality deter- 
mines a chaotic layer. We find that this chaotic layer coincides with the larger one of the 
two unresolved regions in Fig.l. 

2 

B. Case 2, E s2 = The jump of energy can also be rewritten in the same form as in 
Eqs.(24) and (25) . The v c2 is derived from Eq.(18) by assuming that E = E s2 , X = 0. 
The change of phase is derived by evaluating the following integral, 

4> = + lo f ° dx/^2E + 2/{x 2 + a) 1 2, (27) 

J 

where xo is the initial position which is given by Eq.(7). As E — > E s2 , X — > oo and the 
following approximation is applicable, 

e 2 

X ~ -l/2Fcostut , E ~ — - + 2Fcoscut . (28) 
Then we obtain approximately X ~ Eo } E - ■ As E — > E s2 , we approximate Eq.(28) by 

^ ^ V2-E S 2 E s2 — E ^ ^ 

Similar discussions as in case 1 can be done. We obtain the following expression which 
determines a chaotic layer, 

0<£ s2 -£ <(^§M. (30) 

Figure 6 shows that this chaotic layer coincides with the smaller one of the two unresolved 
regions in Fig.l. Actually, in case 2, the change of energy SE s2 for the first collision is large 
enough . This fact makes it possible that these points will fall into irregular region one 
finally That is, in the Poincare section, the two chaotic layers expressed by Eqs. (27) and 
(31) respectively are located in the same connected region. 
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V. ATI SPECTRA AND HOHG IN THE RESCATTERING PROCESS 



Because the amplitude of the field varies slowly with time in the tunneling regime , the 
phase dependence of the ionization rate can be expressed as follows' 8 ! , 

r(0) = lexp(--^— ) a.u.. (31) 
3 3\ecos(p\ 

In terms of the original model in which the Coulomb potential and rescattering are neglected, 
the ATI energy of an ionized electron in linearly polarized laser beams takes the form 

E ATI (4>) = 2U p sin 2 (t). (32) 

Then the energy distribution of the ionized electrons is evaluated by the following expression, 

P(E ATI )dE ATI = r(0)#, P{E ATI ) = r(0)— ^— . (33) 

ati A Ti 

The plots of the distribution of electron energy are shown in Fig. 7, which represents a steeply 
decreasing curve. Some ionization experiments are performed in the tunneling regime to test 
and verify the quasi-static theory. Analysis shows that qualitatively the prediction of the 
two-step theory agrees with the experimental data, and it gives a much smaller fraction of 
the high-energy ionized electrons. M 

In terms of the improved model, because of the rescattering effects, the final energy of 
the detected electrons in a short pulse experiment can be larger than 2U p . The initial phase 
may be divided into four different regions: [1.57,3.32], [3.32,4.05], [4.05,4.52] and [4.52,4.71] 
for uj = 0.07a.w,e = 0.06a.u., and [1.57,3.30], [3.30,3. 70], [3.70,4.69] and [4.69,4.71] for lu = 
0. 04242a. u., e = 0.06a.u.. The trajectories initiating in different regions show quite different 
behaviors (Fig. 8). The electrons starting from the first region drift away from the nucleus 
without a collision. Those electrons originating from the third region must recross the 
original point once before completely escaping. However, the electrons initiating in regions 
2 and 4 can undergo multiple collisions with the nucleus and absorb more photons in this 
rescattering process. As was mentioned above, the trajectories starting from regions 1 and 
region 3 are regular, and from regions 2 and 4 may be irregular. In the situation that the 
initial phase is in region 4, the tunneling barrier is too broad for any significant release to 
occur. Then the effects of chaotic rescattering in region 2 may influence strongly the final 
energy distribution. 

Because of chaotic rescattering the function E AT i is not differentiable almost everywhere 
in those irregular regions . The usual Riemann integral can not be applied to obtain the 
energy distribution, and the Lebesgue integral is introduced to deal with the problem. The 
probability of the kinetic energy of electrons falling into the interval [E, E + AE] is 

P = Les. f E ATI {4>)T{<t>)d<t>, (34) 
J s 

where S denotes a set of initial phase satisfying E < E ATI {<p) < E + AE, AE = U p /10. 
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FIG. 7. Energy spectra for a)F = 0Ma.u.,uj = 0.04242a.u. b)F = 0Ma.u.,uj = 0.07a.u.. 
The solid lines and dotted lines denote the results from the original model and the present model 
respectively. 

In Fig. 7 we plot the energy distribution according to the above formula. It shows clearly 
that rescattering increases the fraction of the electrons with higher energy. This is because 
the electrons which have more probability of staying in the vicinity of the nucleus will 
absorb more photons. This result is also qualitatively in agreement with that of our quantum 
approaches '8]. In numerical simulations, 1000 initial tryjectories are used. We also calculate 
the ATI spectra by using 10000 trajectories . We find that they are agreeable except for a 
little fluctuation in the higher energy range. 

Another consequence of the electron- ion interaction is emission of light. The emission 
can be calculated from the expectation value of the dipole operator 

d{t) =< ip\ex\i/j >, (35) 

where ip is the Floquet state. 

In terms of the symmetry of the Hamiltonian Hq(x,p) = Hq(—x, —p) and the periodicity 
of the external field, Floquets theorem guarantees that the Fourier decomposition of d(t) 
contains neither even nor constant harmonic but odd harmonics dipole component. ' 16 1 

After tunneling the density of states available to electron is large enough so that its 
subsequent evolution should be accurately characterized by classical motion. Then the 
spectral analysis of the classical trajectories can reveal the mechanism responsible for high- 
order harmonic generation in the rescattering process. We consider the power spectrum or 
spectral density of the dipole moment function x(t) in the classical case, 

I(uj) = (l/27r)^m T ->oo(l/2T) < | / ex(t)exp(-iu;t)dt\ 2 >, (36) 

Jo 

where <> indicates an average over field phases <fi. 

Actually analysis of a single trajectory is sufficient to illustrate the underlying 
mechanisms J 17 ' 18 ] Four types of classical trajectory related to rescattering are shown in 
Fig. 8. For regular trajectories (types 1 and 3) the power spectrum consists of only a single 
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peak precisely located at the laser frequency because Rayleigh scattering is dominant in this 
case. Things are much different for trajectories originating from in regions 2 and 4. Figure 
9 shows the dipole moment function and its corresponding power spectrum for a trajectory 
initiated in the irregular regions. In this case , the electron is first pumped by the external 
fields far away from the nucleus for a long time before returning to the ion and emits light. 
Only when the electron is close to the nucleus, can the electron exchange energy efficiently 
with the electromagnetic field. After oscillating around the nucleus for about 10 optical 
cycles, the electron absorbs more photons and is ejected to the continuum. A power spec- 
trum analysis of the dipole moment is made in the time interval [10000, 11000] when the 
electron oscillates near the nucleus. The power spectrum shows that high harmonics up to 
33rd are generated in the process. The magnitude of the peaks decreases monotonically ap- 
proximately with increasing order. The peak corresponding to 33rd harmonic is about three 
orders lower than that of Rayleigh scattering. The even harmonics component appearing in 
the power spectrum can be removed by averaging over an ensemble of classical trajectories. 
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FIG. 8. Typical trajectories initiated in four different regions. 
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FIG. 9. The dipole moment as a function of time and corresponding power spectrum for a 
typical trajectory initiated in irregular regions. 4>q = 3.3250968504, F = 0.06a.u.,uj = 0.07a.u.. 

In the high-intensity regime a new contribution to harmonic emission becomes possible. 
That is , while an electron is near the nucleus it can make a transition back to the ground 
state, emitting a high-energy photon. Obviously the energy of the resulting photon cor- 
responds to the energy (return energy E r ) the electron has when it collides with the core. 
This return energy depends strongly on the field phase at which the electron escapes through 
the potential barrier. In Figs. 10(a) and 11(a) we plot the return energy vs the phase for 
uj = 0. 04242a. u. and uj = 0.021a. -a. respectively. The 1000 initial phases distributed equally 
in the interval [7r/2,37r/2] are used in the calculation. In Figs. 10(b) and 11(b) we show the 
histograms of the return energies weighted by the tunneling rate appropriate to the initial 
conditions for each trajectory. The result is a broad, flat distribution followed by an abrupt 
cutoff. The sharp cutoff in the electron energy occurs at 3.5U P for uj = 0. 04242a. u. and at 
3.3U p for uj = 0.021a. -a. respectively. We think this is the physical origin of the cutoff law 
for high harmonic radiation, t 6 ' 19 ] The smaller the field frequency, the more the tunneling 
condition is satisfied. Therefore the value of the cutoff energy tends to 3.2U P , which has 
been supported by experimental evidence. 
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FIG. 10. a) The energy of returning electron as a function of field phase, b) Energy histogram 
for electron trajectories returning to the nucleus, uj = 0. 04242a. u. and F = 0.06a.u.. 
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FIG. 11. a) The energy of returning electron as a function of field phase, b) Energy histogram 
for electron trajectories returning to the nucleus. u> = 0.021a.w. and F = 0.06a.u.. 



VI. DYNAMICAL STRUCTURE IN POINCARE SECTION 

Irregular behavior of the electrons in the rescattering process is characterized by the 
singularity of the final ATI energy distribution. This singularity is attributed to chaotic 
scattering in system (15). The system studied in previous sections constitutes an excellent 
physical system for the investigation of chaotic scattering. The physical process that we 
have in mind involves the rescattering in the OFI problem. 

To study chaotic scattering, we choose the initial condition as the points in phase space 
well in the asymptotic region of the potential, typically, X = 100 a.u. and P in the range 
[0, —1.5] a.u.. Figure 12 shows the stroboscopic map, (X(t), P(t))\ U!t= o(mod27r), for such an 
ensemble of 200 scattering trajectories. Further the system can flow arbitrarily between the 
bounded and unbounded regimes, so it is a non-compact system. The phase section shows a 
rich and complicated structure. As the system is periodically perturbed, the Poincare section 
is naturally defined as X : cut = 2kn. The corresponding Poincare mapping is expressed as 

M : (X', P') — M o (X, P) = (F(X, P), G(X, P)), (37) 

where (X', P') is the point at ut = ut + 2n of the trajectory initiating in (X, P) at ut . Its 
corresponding tangent mapping will be 

TM : (5X>, 5P>) = Q5X + ^5P, ^5X + ^5P). (38) 

In the following discussions t is set to be zero without losing generality. 
The varitional equations of the Hamiltonian system are 

d(SX) d 2 K d*K 
dt dPdX dPdP ' 

d(SP) d*K d*K 

— = ~dxdx 6X - dxdP 6R (39) 
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Equations (16) and (40) constitute a four-dimensional non-autonomous dynamical sys- 
tem. Tracing the trajectories initiating in (X , P , 1, 0), to — and (X , P , 0, 1), to — by us- 
ing a 4th order Runge-Kutta algorithm in x and p, we can get the solutions (X 1 , P 1 , 5X ± , 5 Pi) 
and (Xi, Pi, SX 2 , SP 2 ) at time ut = 2ir respectively . The Pioncare mapping and its tangent 
mapping can be calculated numerically according to following expressions: 

X l = F(X , P ), Pi = G(X , P ), (40) 
dXi OF 



dX dX 
dPi dG 



sx u 



dX ~ dX 
dXi _ OF 
~dPo~~dPo 
dPi _ dG 
dP ~dP ~ 



= SPl (41) 
= 5X 2 , 

SP 2 . (42) 




-3 1 1 1 1 1 1 1 1 ' 1 
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FIG. 12. Stroboscopic map for 200 trajectories originating asymptotically at X = 100 with P 
in range [0, —1.5]. F = 0.06a.u.,uj = 0.07a.u.. 
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FIG. 13. Phase plane of Poincare mapping M for F = 0.06a.u.,oo = 0.07a.u.. 
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Since the Cartan differential 2-form is conserved under Hamiltonian flow, the map- 
ping M should conserve dX A dP, that is , it is area-preserving. Because of the sym- 
metry of the Hamiltonian system, the 1-periodic fixed points are located on the X-axis. 
The Newton-Raphson algorithm is employed to locate these fixed points precisely. They 
are Pi(- 15.95186600603247,0), P 2 (-15.02, 0), P 3 (-12.25, 0), P 4 (-9.30,0), P 5 (2.80,0), 
Pe(15.87, 0). Pi and P3 are stable elliptic fixed points, others are unstable hyperbolic fixed 
points. Figure 13 shows the phase plane of the mapping M. Around P3 there is a larger 
stable region where KAM curves exist. The stable region around P 2 is too small to be 
distinguished in the plots. 

From the symmetry K{X, — P, —t) = K(X, P, t), the stable manifolds and unstable man- 
ifolds are symmetric about the X-axis. In Fig. 14 we plot the unstable manifolds of these 
four hyperbolic fixed points by using the method suggested in Ref.[21]. Clearly , they inter- 
sect the X-axis. So we conclude that the unstable manifolds and stable manifolds intersect 
each other transversely. The intersection points constitute a hyperbolic invariant set with 
a fractal dimension which has been mentioned in section III. The Smale structure occurs in 
the phase plane. The is the dynamical source of chaotic scattering. Comparing Figs. 12 ,13 
and 14, we find that the unstable manifolds organize the scattering pattern and the phase 
plane structure. 




FIG. 14. The unstable manifolds of period-1 saddles of mapping M for 

F = 0Ma.u.,LO = 0.07a.u.. 



VII. DISCUSSIONS AND CONCLUSIONS 

Recent advances in classical nonlinear dynamics and chaos have had important appli- 
cations in the description of the photoabsorption spectrum of Rydberg atoms in strong 
magnetic fields, ' 22 ^ and the microwave ionization of highly excited hydrogen atoms, t 23 ' These 
discussions extremely enrich the atomic physics. More recently, some pioneer nonlinear dy- 
namicist begin to deal with the topics of atoms in intense fields. Sundaram and Jensen^ 24 ! 
connect the ionization suppression with 'scarred' wave functions. Benvenuto et al} 25 ^ esti- 
mate the parameter regime of atom stabilization according to the border of chaos. Cocke and 
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Reichl^ 26 ! calculate the HONG in the frame of ID model atom and conclude that it implies 
chaos. Richards' 27 ' has studied ID excited hydrogen atoms in short-pulsed high-frequency 
fields and find the signs of incipient chaos. The purpose of this paper is to apply the nonlin- 
ear theory in ATI and HOHG according to a ID description of optical field ionization in the 
regime of low field frequency. We discover a new physical process - chaotic scattering under- 
lying the optical field ionization process. Modern nonlinear theory is applied to analyze the 
chaotic behavior of this nonlinear system. Chaotic layers and self-similar fractal structure 
are found in this model. Otherwise we find that unstable and stable manifolds organize the 
chaotic scattering pattern. 

In another aspect, the above-threshold ionization energy spectra and high-order har- 
monic generation are discussed according to the improved quasi-static model. The rescat- 
tering effects on ATI spectra and high-order harmonic generation are discussed. The result 
that rescattering increases the fraction of the electrons with higher energy is qualitatively 
in agreement with our previous works from a quantum approach. The cutoff law for high 
harmonic radiation is identified by our new model. In particular, We have shown the con- 
nection between regular and irregular trajectories and ATI , HOHG phenomena. We find 
that chaotic behavior plays an important role in these processes. 

In this paper, our discussions are restricted to ID model . A question that can not be 
evaded is: Whether the above discussions are helpful to the understanding of the ionization 
process of a real atom in a three-dimensional (3D) space ? The Key point of this problem is 
whether multiple collisions ( resulting in chaos) still occur in a 3D case. It seems that there 
are only trivial number of multiple return trajectories. Actually, this is a specious argument. 
In our 3D extension of the above model, numerical simulations show that there do exist a 
number of (not trivial number) trajectories which return to the core more than once and 
those multiple collisions trajectories (chaos) play an important role in energy spectra and 
angular distribution of the final hot electrons (results in details will be reported elsewhere) . 
Therefore, we believe that our discussions on ID model is much helpful to the understanding 
of the ionization process of a real atom in a 3D case. 

Notes: After we completed this manuscript, we find that the authors of Ref. [28] have 
made some researches on the effects of multiple returns on double ionization for 3D model, 
and show an enhancement by more than an order of magnitude. 
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CAPTIONS OF FIGURES 

Fig.l Plots of electron energy vs field phase. Dashed lines and solid lines denote the 
results from original model and the present model respectively, a) F — 0.0Qa.u.,u = 
0.07a.u. b)F = 0.06a.u.,uj = 0.04242a.w.. 

Fig. 2 Trapping time vs field phase. F = 0.06a.u.,u = 0. 04242a. u.. 

Fig. 3 Successive magnification of the plot for F = 0.06a.u.,uj = 0.07a.u.. 

Fig. 4 Calculation of fractal dimension in the case F = 0.06a.u.,uj = 0.07a.u. 

Fig. 5 Probability of multiple collisions obtained from making statistics on 4000 tra- 
jectories corresponding to different initial phase (circle). The straight line denotes the 
scaling law mentioned in the text. 

Fig. 6 The chaotic layers predicted theoretically showing coincidence with those unre- 
solved regions. F = 0.0Qa.u.,u = 0.07a.u. 

Fig. 7 Energy spectra for a)F = 0.06a.u.,uj = 0. 04242a. u. b)F = 0.06a. -a., uj = 0.07a.u.. 
The solid lines and dotted lines denote the results from the original model and the 
present model respectively. 

Fig. 8 Typical trajectories initiated in four different regions. 

Fig. 9 The dipole moment as a function of time and corresponding power spectrum for a 
typical trajectory initiated in irregular regions. O — 3.3250968504, F = 0.06a.u.,cu = 
0.07a.w.. 

Fig. 10 a) The energy of returning electron as a function of field phase, b) Energy 
histogram for electron trajectories returning to the nucleus, uj = 0. 04242a. u. and 
F = 0.06a.w.. 

Fig. 11 a) The energy of returning electron as a function of field phase, b) Energy 
histogram for electron trajectories returning to the nucleus, uj = 0.021a.w. and F = 
0.06a.-u.. 

Fig. 12 Stroboscopic map for 200 trajectories originating asymptotically at X = 100 
with P in range [0, —1.5]. F = 0.06a. it., uj = 0.07a. -a.. 

Fig. 13 Phase plane of Poincare mapping M for F = 0.0Qa.u.,uj = 0.07a.u.. 

Fig. 14 The unstable manifolds of period-1 saddles of mapping M for F = 0.06a.w., u = 
0.07a.w.. 
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